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ABSTRACT 

Merging compact binaries are the most viable and best studied candidates for gravitational wave 
(GW) detection by the fully operational network of ground-based observatories. In anticipation of the 
first detections, the expected distribution of GW sources in the local universe is of considerable interest. 
Here we investigate the full phase space distribution of coalescing compact binaries at z = using 
dark matter simulations of structure formation. The fact that these binary systems acquire large 
barycentric velocities at birth ("kicks") results in merger site distributions that are more diffusely 
distributed with respect to their putative hosts, with mergers occurring out to distances of a few Mpc 
from the host halo. Redshift estimates based solely on the nearest galaxy in projection can, as a result, 
be inaccurate. On the other hand, large offsets from the host galaxy could aid the detection of faint 
optical counterparts and should be considered when designing strategies for follow-up observations. 
The degree of isotropy in the projected sky distributions of GW sources is found to be augmented 
with increasing kick velocity and to be severely enhanced if progenitor systems possess large kicks 
as inferred from the known population of pulsars and double compact binaries. Even in the absence 
of observed electromagnetic counterparts, the differences in sky distributions of binaries produced by 
disparate kick- velocity models could be discerned by GW observatories, within the expected accuracies 
and detection rates of advanced LIGO-in particular with the addition of more interferometers. 
Subject headings: gravitational waves — stars: neutron — binaries: general 



1. INTRODUCTION 

The merger of double compact objects represents the 
first identified and most predicta ble source of gravita- 
tional wave (GW) radiation (e.g. iPhirm ev 1991). Only 
recently have the first GW observatories come online, 
and the first detection events are expected in the next few 
years. Over the past three decades the merger rate within 
the local universe has been thoroughly examined (see e.g. 
lAbadie et all 120101 : iMandel fc O'Shaughnessvl 120101 for 
recent reviews). 

The merger rates are expected to be dominated 
by mergers of neutron-star binaries, with (-ft) ~ 
lMpc _3 Myr _1 . However, these rates are signifi- 
cantly uncertain, since they come either from extrapo- 
lations from the small observed sample of Galactic bi- 
nary pulsars whose luminosity distribution is not well 
constrained, or from population-synthesis models that 
have many ill-determined paramet ers such as cqmmon - 
envelope efficiencies. In particular, [Abadie et al.l (|2010T ) 
estimate the confidence bounds on the neutron-star bi- 
nary merger rates as 5ft « 0.01 — 10Mpc~ 3 Myr~ . 
The horizon distance^ for the initial and advanced 
LIGO/ Virgo detector networks are estimated as T> ~ 30 
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and ~ 400 Mpc, respectively, based on the distance at 
which a single detector could detect gravitational waves 
from a neutron-sta r bina ry at a signal-to-noise ratio 
of 8. lAbadie et al.l (j2010l ) estimate that the advanced 
LIGO/Virgo network could plausibly detect between 0.4 
and 400 neutron-star binaries per year, with a likely rate 
of approximately 40 detections per year. 

The prospects for detection and characterization of 
GW sources are thus sensitive to the distribution of com- 
pact binaries in the local universe (i.e. at distances 
< T>). The fact that thes e systems must have large sys- 
temic velocities at birth (jBrandt fe Pods iadlowskii [i995l ; 
iFrver fc Kalogeral I1997D implies that by the time they 
merge, after approximately a Hubble time, they will be 
far from their birth sites. The locations of merging sites 
depends critically on the binary's natal kick velocity and 
the temporal evolution of the gravitational potential of 
the host halo as well as that of its nearby neighbors 
(|Zemp et al.ll2009l ). 

In this Letter, we study the evolving distribution of 
compact binary systems from formation until coalescence 
at z = using cosmological simulations of structure for- 
mation. This allows us to examine the full radial and 
angular distributions of merging compact binaries in the 
local universe. The organization is as follows. In ^ 
we describe the numerical methods and initial setup and 
the criteria used to select a local-like universe. The dis- 
tributions of compact binaries at z — are presented in 
|}3] for three different kick velocity scenarios, and in i|4] 
we examine the ability of GW observatories to discern 
between them experimentally. Finally, discusses the 
implications of our findings. 
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2. METHODS AND INITIAL MODEL 
2.1. Simulation 

The focus of this work is to understand the distribu- 
tion of compact binaries in the local universe using cos- 
mological simulations. To this end, we have performed a 
dark matter only cosmological structure formation sim- 
ulation following the numerical procedure outlined in 
iZemp et all §009). A comoving 80 Mpc periodic box 
is initialized at redshift z = 22 A (161 Myr) and use s 
WMAP3 cosmological parameters (Spergel e t al . 1 120071 ) . 
The initial conditio ns are evolved using the parallel tree 
code PKDGRAV2 (jStadell [200l until z = 1.60 (4.24 
Gyr). At this time, we populate each halo with mass 
greater than 2.15 x lO 11 ^/© (of which there are 2461 in 
the simulation) with 2000 massless tracers. 

Each tracer particle is meant to represent a compact 
binary system, which, on av erage, forms around the p eak 
of the star formation epoch ()Madau et al.lll996lll998f ) . In 
general, the local merger rate is given by the convolution 
of the star formation rate with the probability distribu- 
tion P(t) of the merging time delays r. Compact binaries 
formed at the peak of the star formation history, merg- 
ing after delays consistent with the orbital separations of 
known relativistic binary pulsars (|0'Shaughnessv et al.l 
2008), dominate the local merger rat<S 

Tracers are injected into the center of their halo, 
with an isotropic Maxwell-Boltzmann velocity distribu- 
tion with mean speeds v — 360, 180, and 90 km s _1 and 
dispersions a = 150, 75, and 37.5 km s _1 (hereafter 
denoted as models M 360 , Mi 80 and M go ). This is con- 
sistent with the magnitude of the natal kicks required to 
explain the observed parameters of binary NS systems 
([Brandt fc PodsiadlowsklH995l : iFrver fc Kalogeraill997H 
— only when the center of mass kicks have magnitudes 
exceeding 200 km s _1 can the progenitor orbits be suf- 
ficiently wide to accommodate evolved helium stars and 
still produce the small separations measured in these sys- 
tems. 

The contribution of individual tracers to the over- 
all population is weighted linearly with their progenitor 
halo's mass (at z — 1.6) in all of our calculations. 

Finally, the cosmological box together with the tracer 
particle populations are evolved until redshift z = (13.8 
Gyr). This results in diverse predictions of compact bi- 
nary demographics at z = in the case of an isotropic 
kick velocity distribution whose properties do not vary 
with the initial binary separation. Merger times in pop- 
ulation synthesis models are found to be relatively in- 
sensi tive to the initial kick velocity (e.g., iBloom et al.l 
119991) . This not only justifies our assumption but, when 
taken together with the pr ogenitors' long time delays 
(jO'Shaughnessv et al.ll20lOT ). also reinforces the validity 
of a single injection time. 

2.2. Local-Like Universe Selection 

Once the tracers and DM are evolved to z = 0, a local- 
like universe is select ed based on criteria adapted from 
IHoffman etall (pOOl . The local-like universe is charac- 
terized here by the following: 

6 For P(t) oc 1/r this early-assem bled population c ould increase 
the local event rate by at least ~ 3 (Guctta & Piran 2005). 



i. There are two dark matter halos, representing the 
Milky- Way and Andromeda pair, with maximum 
circular velocities V c G [125,270] km s _1 . 

ii. These halos are separated by d < 1.4 h^ 1 Mpc, and 
approaching each-other (i.e. d < 0.0 km s _1 ). 

iii. There is a Virgo-like halo at a distance d G [5, 12] 
h^Mpc, and V c G [500, 1500] km s" 1 

iv. No halos with comparable or higher maximum cir- 
cular velocities than either of the pair exist within 3 
/i _1 Mpc, and no other Virgo-like halos exist within 
12 h^Mpc. 

The first three constraints resulted in three local-like 
groups. Inclusion of the fourth criterion resulted in a sin- 
gle, optimal environment for our analysis. These criteria 
ensure that the evolutionary environment is similar to 
that of the actual Milky Way and local group galaxies. 

3. THE LOCAL DISTRIBUTION OF COMPACT BINARIES 

We now examine the local, three-dimensional distribu- 
tion of merging compact binaries, which are character- 
ized here by the massless tracer particle population at 
z = centered on the Milky Way-like galaxy as defined 
in £12.21 Figure [T] shows the radial distribution of trac- 
ers and dark matter within our local-like universe. In 
models Migo and Mgo, tracer particles closely follow the 
dark matter central-density pe aks, just like the galaxie s 
themselves in CDM cosmology (Blu menthal et aLlll984l) . 

As the natal, barycentric kick-velocity becomes com- 
parable to the escape velocity of the progenitor halos, 
an increasing fraction of tracers escape. These unbound 
tracers pollute the intergalactic region, thus forming a 
particle background which closely follows the overall dark 
matter distribution, as seen for model M360 in Figure [TJ 
In super-galactic regions with more continuous gravita- 
tional potentials (i.e. more densely and uniformly pop- 
ulated halo environments, such as clusters), the tracer 
background becomes more heavily populated. The rel- 
ative isolation of the Milky Way- Andromeda pair con- 
tributes to the enduring presence of strong native tracer 
peaks located at each halo center. As a result, the ex- 
tended background distributions of tracers-centered on 
the pair-are only apparent in the highest kick velocity 
model. 

As seen in Figure [3J the number of tracers within a 
sphere encompassing the Milky Way and Andromeda ha- 
los is noticeably depleted at higher kick-velocities. At 
these velocities, the host halos are unable to effectively 
retain most of their tracers. As the sphere's volume ap- 
proaches the Virgo cluster, the number of (sub)halos be- 
comes so large that the mean separation between cen- 
tral peak densities decreases below the characteristic size 
of the background tracer population. As a result, the 
variation with kick-velocity in the tracer distributions is 
drowned-out. It should be noted that the e xpected event 
rate in such a small volume is negligible (lAbadie et al.l 
2010), and thus the effects of varying kick-velocity will 
be indiscernible in the integrated merger-rate of compact 
objects within LIGO/Virgo detection horizons. 

Although the integrated tracer distribution is insen- 
sitive to the model, the angular distribution of tracers 
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depends strongly on the binary's kick velocity. This is 
evident in Figure [3j which plots sky maps of tracers and 
dark matter within a given volume (resolved to 4 square 
degree pixels). As expected, high velocity kicks lead to 
more pronounced isotropies when compared to the low 
velocity kick scenarios. At 10 Mpc, ~40% of the M360 
weighted tracers lie in pixels outside those of M90; this 
fraction falls to 15% and 10% for 40 and 80 Mpc respec- 
tively. This trend results from the increasing isotropy of 
dark matter in projection at progressively larger scales. 

For large velocities, the distribution of GW sources 
forms a sky continuum (Figure [3]) rather than well- 
isolated substructures — complicating host galaxy iden- 
tification and thus redshift determination. On the other 
hand, the extension of the particle tracer distributions, 
which grows with increasing kick velocity, could aid the 
detection of photonic counterparts, especially at optical 
wavelengths. This is because at large kick velocities the 
majority of the mergers will take place well outside the 
host galaxy's half-light radius. 



4. PREDICTIONS FOR GRAVITATION WAVE 
OBSERVATIONS 

The number of detections required for GW observato- 
ries to be able to reconstruct the kick-velocity distribu- 
tion is examined here. Timing triangulation from rela- 
tive GW phase shifttQ between widely separa ted detec- 
tors i s the primary source of sky localization (jFairhurstl 
12009( 1 , and Fisher matrix or Markov Chain Monte Carlo 
techniques can be used to compute error estimates 
(|van der Sluvs et all 120081 ). The intricacies of parame- 
ter determination and error estimation can be extensive, 
as correlations between waveform parameters mean that 
some parameters (such as dist ance and inclination) are 
partially degenerate (see, e.g. iCutler fc Flanagan! fl994L 
and references therein). Typically, with a three-detector 
network, the distance to the GW event could be deter- 
mined to within a ^20-50% uncertainty, and the angular 
location of the event to ~ 5 — 50 square degrees, de- 
pendi ng on source location, masses, and signal-t o-noise 
ratio (jFairhurstl [200l Ivan der Sluvs et all I2008I ). The 
large uncertainty in the distance determination can be 
understood when considering its dependence on the sig- 
nal amplitude, which is much more uncertain than the 
phase-space information. 

To estimate the number of events required to dis- 
tinguish between different kick velocity models, we 
apply a Baye sian approach similar to that used by 
Mandcl (2010) to approximate the efficacy of popula- 
tion reconstruction from GW signals. Data sets Dj G 
-D36O; -D18O arL d Dgo are drawn from each model Mi £ 
-^360,-^180 and M90 respectively. Each data set con- 
tains n independently drawn data points (i.e. tracers), 
characterized by 3 position coordinates; i.e. Di{n) = 
[xi t i(r,a,S),Xi t 2(f,Oi,S), . . .Xi tn (r,a,5)}. The probabil- 
ity of a tracer being selected for a given data set is lin- 
early proportional to the halo mass of the progenitor, in 
accordance with the weighting scheme of §2.11 The prob- 
ability that a particular model i fits a data set j can be 



7 Fo r a review of GW emission from compact binaries see lHughesI 
p009h . 



rewritten using Bayes' formula: 

PiD^Mi) ■ P{Mi) 



P(M t \D 3 (n)) 



P(D 3 (n)) 



(1) 



Throughout our analysis we assume flat priors 
[P(Mi) = P(Mj)], and equivalent evidence [P(A) = 
P(Dj)\. A comparison between models then yields: 



P(Mi\Dj(n)) = P(A(n)|Mj) 
P(M,|AH) P(AH|Mj) 



n 

fe=i 



P(x l .k\M i ) 
Pix^Mj)' 



(2) 



where the probability of a particular data point given a 
specific model, P(xi k\Mj), is described by the convolu- 
tion of the point spread function (PSF-5) of the detector 
with the probability distribution function of the model 
in question. That is 



1 

P(xi, k \Mj) =^2S(x i>k \pmel l ) • P(pixel;|M,), 
1=1 



(3) 



where the sum is being performed on each pixel (pixel;) 
for all q pixels. 

The PSFs are assumed here to be gaussian in each 
coordinate direction, characterized by standard devi- 
ations in distance, right ascension, and declination: 
tfhish = [5%, 1°,1°], a mcd = [30%,2°,2°], a low = 
[50%, 4°, 4°]. These reflect different assumptions for 
the high, medium, and low accuracy of positional re- 
construction for gravitational-wave detections. The ex- 
act parameter-estimation accuracy is difficult to predict, 
since it will depend both on the details of the detector 
network (e.g., the relative sensitivity of detectors and 
their calibration accuracy) and on the specifics of indi- 
vidual events (their signal-to-noise ratio, and the sky lo- 
cation and orientation of the binary). Therefore, these 
three assumptions should be considered only as possible 
predictions for typical accuracies. Thus, low accuracies 
may be typical for events detected with a three-detector 
LlGO/Virgo network at the threshold of detectability. 
Meanwhile, the addition of a fourth interferometer, such 
as a possible AIGO detector in Australia or LGCT in 
Japan, could significantly enhance the sky localization 
accuracy and moderately improve distance sensitivity 
(Fairhurst et al. 2010), making medium- accuracy mea- 
surements typical and high-accuracy measurements pos- 
sible. 

In these calculations, we compare hypothetical GW ob- 
servations with models of compact binary distributions. 
This comparison is being made assuming that the local 
dark matter distribution is perfectly known. In reality, 
this is not the case; and the results presented here are 
thus optimistic. In the future, the comparison between 
model and observation should be refined to include the 
local distribution of light (e.g. galaxies) rather than dark 
matter halos. 

Table Q] summarizes the ability of GW observatories to 
discern the kick velocity distribution of the merging bi- 
naries from the reconstructed angular positions and dis- 
tances (assumed to be determined without a galaxy host 
association). Two sample volumes are considered: 40 
and 80 Mpc. This is done in order to understand the 
sensitivity of our results to the uncertainty in physical 
separation which, for a fixed angular resolution, varies 
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with distance. We find that ~ 50 events are required to 
distinguish between the lowest and highest kick velocity 
scenarios for moderate detector accuracies, irrespective 
of which sample volume is examined. For low detector ac- 
curacies, 50 — 350 detections are necessar 50. Thus, a dis- 
tinction between the two extreme models is possible once 
advanced detectors come online, with an expected event 
rate of ~ 40 per year for detect ions of binary neutron 
star mergers (jAbadie et alJl2010[) Fl Meanwhile, distin- 
guishing between the two low-kick scenarios is very dif- 
ficult, if not impossible, until the era of third-generation 
detectors. The addition of a fourth GW detector to the 
LIGO /Virgo network would significantly improve source 
localization, and thereby the accuracy with which event 
distributions could be distinguished. 

Assuming a LIGO/Virgo horizon of ~ 400 Mpc, only 
~ 10% of all detected mergers would take place within 
80 Mpc. With a constant angular resolution, the uncer- 
tainty in physical position is proportional to the event's 
distance, suggesting that using events at greater dis- 
tances leads to a degradation in the ability to distinguish 
between kick-velocity models. Although we find no clear 
increase in the number of required detections between the 
40 and 80 Mpc samples, further investigation is required 
to assess the effects of a larger sample volume. 

5. SUMMARY 

In this Letter, we use dark matter cosmological simula- 
tions to examine the full three-dimensional distribution 
of coalescing compact binaries in the local universe un- 
der the following assumptions. First, we assume a sin- 
gle epoch of star formation and a simple star formation 
recipe; that is, the contribution of a particular halo to 
the total star formation is directly proportional to its 
dark matter mass. Although a more realistic treatment 
of star formation should be considered, we do not ex- 
pect that our qualitative results will change significantly. 
Second, we assume an isotropic natal kick velocity dis- 
tribution, whose properties are invariant of initial binary 
separation. Under this assumption, the merging time is 
independent of the kick velocity. This is found to be a 
reasonable approximation in binary population synthe- 
sis models, which helps justify our single epoch of tracer 
injection. Third, our comparisons between kick veloc- 
ity models in <j4] assume a perfect knowledge of the local 
dark matter distribution, when in actuality this distribu- 
tion would have to be deduced from the observable, lo- 
cal universe. Finally, due to computational constraints, 
only an 80 Mpc region of the expected 400 Mpc hori- 
zon of advanced LIGO/Virgo has been modeled. Despite 
the increased uncertainty in the true-distance offset be- 
tween host and merger at such distances, the difference 
between our 40 and 80 Mpc results (Table [T]) suggests 
that our methods could remain effective in deducing the 
kick velocity distribution with a reasonable number of 
detections. Keeping these assumptions in mind, it is still 
evident that the use of static, non-evolving potentials for 
individual hosts at the time of binary formation severely 

8 It is important to note that the number of detections required 
is highly sensitive to the model from which the data is drawn, not 
simply on which models are being contrasted. 

9 The event rate estimates have significant uncertainties, and 
range from pessimistic estimates of ~0.4 events per year to opti- 
mistic estimates of ~400 events per year (Abadic et al. 2010). 



overestimates the ret ention of all but the lowest barycen- 
tric velocity systems ([Fryer et al.lll999t iBelczvhski et al. 
2000tlRosswog et alj2003tlBloom et al.ll999H B ulik et al. 
1999 : iPortegies Zwart fc Yungelsonlll998h . 

Static calculations predict that the distribution of 
gravitational wave sources in the sky should closely trace 
the distribution of galaxies. An accurate inclusion of 
evolving host halo potentials in cosm ological simulation s 
have shown this to be inaccurate (|Zemp et al.l [2009) . 
In fact, we show that not only do the distributions of 
merging compact binaries extend well beyond their birth 
halo, but variations in kick velocity lead to marked dif- 
ferences in their sky distributions. The repercussions of 
this result are twofold. On one hand, we find that the 
variation in the projected distribution of double com- 
pact objects with different natal kick-velocities should 
be distinguishable with the expected accuracies of GW 
observatories. In principle, this will allow important 
information on the formation and evolution of the bi- 
nary progenitor to be deciphered from the distribution 
of GW detections alone. On the other hand, the fact 
that the distribution of merging binaries does not ac- 
curately trace the locations of their birth halos com- 
plicates redshift determination. Having said this, the 
presence of a binary distribution extending well beyond 
the half-light radius of their hosts suggests that associ- 
ating optical counterparts to GW events could be eas- 
ier as they are less likely to be drowned out by their 
host galaxy's light. This is particularly important as the 
optical counterparts a r e predicted to be relatively dim 
dLi fc Paczvnskil 119981: iRosswog fc Rarnirez-Rur3 l2002t 
lKulkarnil l2005; Mc tzger et al.ll2010D . 

Gravitational waves offer the possibility of casting 
proverbial light on otherwise invisible phenomena; they 
will - by their very nature - tell us about events where 
large quantitites of mass move in such small regions that 
they are utterly opaque and forev er hidden from direct 
electr omagnetic probing (see, e.g. iLee fe Ramirez-Ruizl 
l2007f h A time-integrated luminosity of the order of a 
fraction of a solar rest mass is predicted from merging 
compact binaries. Ground-based facilities, like LIGO, 
GEO600 and Virgo, will be searching for these stellar- 
remnant mergers in the local universe. The distribution 
of merger sites is thus of considerable importance to GW 
observatories. The proposed use of galaxy catalogs as 
priors when passing triggers from possible GW detec- 
tions to point telescopes for electromagnetic follow-ups 
will need to account for the possibility of mergers away 
from the observed galaxies. Using cosmological simula- 
tions of structure formation, the local sky distributions 
are found to vary with the kick velocity distributions of 
the progenitor systems, allowing a determination of the 
cosmography of massive binary stars. Despite the fact 
that individual detections lack the positional accuracy 
of electromagnetic observations, it may be possible to 
strengthen the case for (or against) high natal kick ve- 
locities based solely on GW observations. The addition 
of more gravitational- wave detectors to the LIGO /Virgo 
network will greatly improve our ability to distinguish 
between models with different kick velocity distributions 
by improving the positional reconstruction of individual 
events. 
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A^360(-D360) 


Mi 80 (Di8o) 


M go {Dgo) 


Dist 


PSF Accuracy 


VS. Mi 80 


VS. Mgo 


vs. M 360 


vs. M 90 


vs. M 360 


VS. Mi 80 


< 80Mpc 


High 


22 


16 


26 


> 1000 


22 


282 


Med 


73 


39 


35 


> 1000 


31 


384 


Low 


> 1000 


349 


52 


> 1000 


50 


881 


< 40Mpc 


High 


27 


17 


34 


> 1000 


23 


> 1000 


Med 


78 


46 


40 


> 1000 


37 


> 1000 


Low 


146 


137 


56 


> 1000 


56 


> 1000 



TABLE 1 

Number of detections required to achieve 99% confidence in the correct model for 90% of data sets. These results are compared 

between two different sample radii, and three different detector accuracies characterized by standard deviations (in distance, right 
ascension, declination) of: a high = {5%, 1°, 1°}, <? m ed = {30%, 2°, 2°}, ai ow = {50%, 4°, 4°}. Entries marked '> 1000' failed to reach 
the desired confidence in the required number of data sets within the 1000 data points used. 
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Fig. 1. — Tracer vs. dark matter distribution in a local-like universe as a function of barycentric kick- velocity. Integrated particle mass 
in uniform radial-width shells is plotted versus distance from a solar-equivalent offset from the Milky Way center. The vertical axes are 
plotted in arbitrary units of number per unit length, with tracers normalized with respect to the total tracer population as described in 
£|2,ll As the kick-velocity increases from 90 km s _1 (top panel) to 360 km s -1 (bottom panel), a noticeable portion of tracers becomes 
dclocalized, forming a background and mixing populations. 
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Fig. 2. — Cumulative distribution of tracers and dark matter within a given distance as a function of kick-velocity. The vertical axes 
are plotted in arbitrary units of number per unit volume, with tracers normalized self-consistently. Although the number of tracers in 
the central halo is noticeably lower for the highest kick-velocity model, the difference is negligible once the volume reaches the Virgo-like 
cluster, where the background distribution of tracers outweighs changes in local distributions. The number of merging binaries is assumed 
to be proportional to the mass of the host halo. 
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Fig. 3. — Sky maps of dark matter (first column) and tracers with highest and lowest kick velocity scenarios (second and third columns, 
respectively) as a function of distance. Figures make use of Hammer projections with 2°x 2° bins. Densities arc plotted in units of column 
density, scaled to the maximum densities of the dark matter and normalized tracer distributions independently. Pixels with no tracers 
or dark matter are white, corresponding to densities less than the resolution of the simulation. Although tracer peak densities remain 
relatively unchanged, a tracer-background forms as the kick velocity approaches the escape velocity. At 90 km s — 1 the tracers follow only 
the dark matter overdensities, as does the light-distrubtion. The distributions approach isotropy slower than the dark matter distribution. 
Differences in distribution are clearly apparent. Note the logarithmic color scale. 



